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1. Introduction 

The evolution of pionium, the hydrogen-like atom formed by a 7r + 7r~ pairs, in a material 
target has been thoroughly studied in the recent years [TJ El HI 03 M Hj due to its 
crucial implications in the DIRAC-PS 212 experiment jHj- This experiment is devoted 
to measure the lifetime of pionium, intimately linked to the strong interaction scattering 
lengths as we will see in section testing the predictions of the Chiral Perturbation 
Theory on these magnitudes. 

The transport of pionium in mater has been always treated using a 
classic probabilistic picture neglecting the quantum mechanics interference between 
degenerated states with the same energy. In the case of hydrogen-like atoms this is 
of particular importance since the accidental degeneracy of the hamiltonian increases 
the amount of states among which the interference can be significative. 

In JO] the density matrix formalism has been used to propose a new set of equations 
for the pionium evolution accounting for the interference effects. In this work we have 
solved these equations and analyzed the consequences for the framework of DIRAC 
experiment. 

2. The Problem of Pionium in Matter 

Due to the short lifetime of the pion, pionium, the hydrogen-like 7t + tt~ atom can not 
be produced at rest in the laboratory frame. However, pionium can be originated in 
collisions of high-energy projectiles with a fixed target. The production cross section is 
given in [TU] : 



where the rightmost term accounts for the production of rc + and tt pairs at equal 
momenta (p = q)- 



eigenstate quantum numbers, n^, and m i} of the hydrogen-like hamiltonian. For 
simplicity, in this work we have chosen to work with monochromatic atoms of 4.6 GeV/c, 
the mean value of laboratory momentum of pionium in DIRAC, moving in the z axis 
direction. The effect of using the experimental pionium laboratory momentum spectrum 
is small as shown in [7] . The yield of a particular state is proportional to its wave function 
squared at the origin. It has been shown ^T] that the effect of the strong interaction 
between the two pions of the atom significantly modifies |"0i(O)| in comparison to the pure 
Coulomb wave function. However, the ratio between the production rate in different 
states has been demonstrated to be kept as for the Coulomb wave functions [T2|. Thus, 
considering that the Coulomb functions obey 




(1) 



The state of pionium is defined by the center of mass momentum P and the 




(2) 
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we see that only S states are created following the 1/nf law. 

The atom moves in a fixed thickness target disposed in the Oz axis and considered 
infinite in the transverse (x, y) coordinates. The target is made of a chemically pure 
material like Nickel, Platinum or Titanium. Our goal is to know the population 
probability of every bound state as a function of the position in the target, z, and 
from this extract other results as the break-up probability. Usually a classical approach 
is used to solve this problem [21 EJ- It consists of considering the total a\ ot and transition 
between two discrete states cry cross sections for a pionium-target atom scattering and 
apply the probabilistic evolution equation: 

= —TiPiizi-no^c^Piiz), (3) 

where P%{z) is the classical probability for the atom to be in the i state, /3j = 16.48 
the Lorentz center of mass to laboratory factor for P = 4.6 GeV/c, no is the number of 
target atoms per unit of volume, and Cij are the transition coefficients. 

The value of n is a function of the density of the target, p, the Avogadro number, 
N , and the atomic mass of the target atoms, A: 

P N ° (A\ 
no = (4) 

while the transition coefficients depend on the pionium-target atom cross sections as: 

(k,i = 5 i: ial ot - a i)h (5) 

The pionium decay is strongly dominated (BR > 99% ^U]) by the 7r + 7r~ — > 7r°7r° 
reaction. Taking this into account, the width of the i state is proportional to the isospin 
and isospin 2 pion-pion scattering lengths difference [14] : 



167r ^Ml - - \M W 
r, = -g w K - oo) (1 + *r 



W C) (0) , (6) 



2 



where M v and M n o are the masses of the charged and the neutral pion and 5r = 0.058 
the Next to Leading Order correction that includes the effect of the strong interaction 
between the two pions. Of course, the width of a state holds Ti = where Xj is the 
corresponding lifetime of the state. Due to (j2J) we can see that pionium only decays 
from S states and the lifetime of any S state is related to the lifetime of the ground 
state: 

r„oo = n 3 r. (7) 

The lifetime of pionium is hence the only parameter to be inputed in the evolution 
equation and can be related to any of its outputs. In particular we will link it to the 
break-up probability. The experimental result of DIRAC will be used to test with 5% 
accuracy the accurate Chiral Perturbation Theory prediction of a[J — cJq = 0.265 ±0.004 
which leads to the lifetime value of r = (2.9 ±0.1) • 10~ 15 s |T5j . 
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3. The Density Matrix Evolution Equation 

Equation (jBJ) has been accurately solved obtaining the eigenvalues and the 
eigenvectors [2] and also with Monte Carlo [7] for the bound states with n < 8, which is 
enough to precisely calculate the break-up probability as explained in [Zj. However, the 
work of Voskresenskaya jH] demonstrates that the use of the classic probabilistic picture 
might be inaccurate. This is because (j3J) neglects the quantum interference between the 
pionium states during their passage through the target. 

A more precise description of the system dynamics is given in terms of the density 
matrix p^. The evolution equation in this formalism is given by 



dpik _ J_ 
dz /?7 



lmPlm{z) (8) 



l,m 

where e k indicates the bounding energy of the k state and flik,i m stands for the transition 
coefficients matrix. This equation reduces to (jHJ), identifying pa(z) = Pi(z), if the Qik,i m 
crossed terms obeying i ^ k or / ^ m were zero. 

The goal of this work is to solve this equation and determine how it corrects © for 
the particular conditions of DIRAC experiment, namely for the result of the break-up 
probability. 

4. The Matrix Elements 

To calculate the matrix elements c^j and flik,im we have applied the coherent pure 
electrostatic first Born approximation approach. Even though it is known that 
relativistic and multiphoton exchange must be accounted to achieve the precision of 
1% [Zj our goal was to check wether quantum interference is a relevant effect. For this 
we will show that pure electrostatic first Born approximation is enough. 

The expression for the pionium-target cross sections in the electrostatic first Born 
approximation, used in the classical picture, was obtained by S. Mrowczyhski time 
ago 0: 



tot 



J 2 J \U{q)\ 2 [l-Ft{q)]d\ (9) 



^/'H^(fHH) 



2 



d 2 q } (10) 



where q is the transferred momentum between the target and the pionic atoms. The 
cross section does only depend on the two transverse coordinates of the momentum due 
to the symmetry of the collision with respect to the scattering axis. We have chosen the 
Fourier transform of the target atom potential U(q) to be the Moliere parameterization 
for the solution of the Thomas- Fermi equation [TH] : 

rrt \ a ry ( 0.35 0.55 0.10 \ 0.3Z 1 / 3 , . 

U(q) = AnZa — k + — ^ + — k q = , (11) 

W \q 2 + qi q 2 + 16q% q 2 + 400^ J W 0.885a 
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being = 0.529 x 10 28 cm the Bohr radius of Hydrogen, a the fine structure constant 
and Z the atomic number of the target atoms. The F-(q) are the pionium form factors: 



^(r)e^*(r)df, 



(12) 



calculated in |5j and [3]. In this work we shall use the code of [T7] based on the result 

The equivalent of © for the fli k ,im elements in the density matrix formalism is 
given by: 



n 



ik,lm 



- n {2) 

ik,lm ik,lmi 



(13) 



where: 



n 



(i) 

ik,lm 



>k,m 



2p 2 



+ 



2jP 



\U(q)\ 2 [25 hl -Fl(q)-Fl(-q)] d 2 q + 
\U(q)\ 2 [25 Km - F^q) - i^(-g)] d\ 



(14) 



plays the role of the total cross section, while 



n 



(2) 

ikJm 



1 



IP 



Fill 



x 



(15) 



would be the analogue of the transition cross section. In fact fi^ m becomes the total 



(2) 

cross section if i = k = I = m and lm the transition cross section if i — k and I = m. 

Equations $H| and p5)l are our main tool for the numerical calculations and their 
development from the original formulas of 9] can be followed in Appendix A 



4-1- Selection Rules and Transition Elements Examples 

As pointed out in 0, and due to the properties of the form factors under the parity 
transformation, the flik,im coefficients are different from zero only if: 

mi - m k - mi + m m = , U - l k - k + l m = 2s (16) 

where we should remember that rrii(k,i,m) an d k(k,i,m) are the magnetic and orbital 
quantum numbers of the states \i(k, I, raj). The index s is an arbitrary integral number. 

For the election of the Oz the quantization axis the transitions between 

states of different ^-parity are strongly suppressed This means that only states with 
even I — m will be populated since pionium is produced in S states only. This, together 
with (|16|) means that: 

Pik{z) 7^ if rrii = m k , h = h + 2s. (17) 
This rule could be broken by the complex coefficient in (JH|) : 
.{ £ k ~ £i) 
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Figure 1. The plots shows the solution of ijSJ for P|ioo)<20o| ; oscillating with high 
frequency around 0, and P|600)(700|i constantly over in a much larger range. 



which produces an oscillatory term in the solutions. However, for the ground and lowest 
excited states the condition: 

kfc — £i\ 



n \&ik,ik\ ^ 



Pi 



holds and the Pik(z) solution oscillates many times in a small interval, compared to the 
electromagnetic transition range (given by n \£lik,ik\) an d can De considered to average 
as zero: 

Pik(z) ~ 0. 

There is an exception if the i and k states belong to the same shell since the energy of 
the hydrogen-like system does only depend on the principal quantum number. In this 
case Ek — Si = 0. Hence, for the low energy states, we can complete relation (|T7|) as: 



p ik (z) ^ if Si = Si 



71; 



m k 



h = l k + 2s. 



;is) 



However, if the principal quantum numbers of the i and k states hold > 6, then 

£fc — £ i\ 



n \&ik, 



ik 



Pi 



and the solution for pik{z) is not zero even though i and k are not states from the same 
shell. 

In figure ^ we can see that whereas pik{z) oscillates more than six times in 0.1 pm 
for p|ioo)(2oo| h d° es n °t oscillate at all for P|6oo){70o| i n a wide range. 

As an example of the matrix elements we consider the subspace formed by the 1 2 1 1 ) , 
1 300) and 1 320) states. The Q matrix restricted to this subspace is shown in tabled 
We can see that at least for the 1 320) (300 1 mixed state the matrix elements are of the 
same order of magnitude as for the same shell pure states. 



5. Solving the System 



We have numerically solved the differential equation systems and (jSJ) using the 
Runge-Kutta method [IS) . Finding the eigenvalues, as in [2j, would be too lengthy due 
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Table 1. ft matrix elements in the 1 2 1 1) , 1 300} , 1 320} subspace. Units are 10 20 barn. 



^\i){k\,\l)(m\ 


|211)(211| 


|300)(300| 


|320)(300| 


|320)(320| 


|211)(211| 


-4.66 


0.044 


-0.083 


0.234 


|300)(300| 


0.044 


-18.2 


2.41 


0. 


|320)(300| 


-0.083 


2.41 


-13.7 


2.41 


|320)(320| 


0.234 


0. 


2.41 


-9.10 



to the size of the density matrix system. The Monte-Carlo method of [Hj cannot be 
applied to the density matrix formalism since the system (jHJ) does not obey: 

^ik,ik ^ ^ \^ik,lm\- (-^9) 
Im^ik 

We have considered a Nickel 95 pm target and a monochromatic 4.6 GeV/c atom 
sample. The lifetime of the ground state of pionium was supposed to be T" 1 = 2.9-10~ 15 s 
according to the theoretical prediction jT^j. The initial conditions are given by: 

Pi(0) = P«(0) = nr 3 /C(3) if Z< = 0, 

p.(0) = ^(0) = otherwise, (20) 

here £(3) = J2 n ~ 3 ~ 1-202. The system has been restricted to the bound states with 
n < 7. This means 84 mixed states and 353088 Q matrix elements different from zero. 
Cutting the number of considered states does only slightly affect the solution of the last 
two cores taken into account (in this case states with n = 6 and n = 7) as shown in 

To achieve a very good accuracy in the final results we have considered a sequence 
of step lengths in the numerical integration of the system: 

h = 2- HT 3 , 1 ■ HT 3 , 0.5 ■ HT 3 , 0.25 ■ 10" 3 , 0.1 • 10~ 3 [/H 

and made a polynomial extrapolation to the limit h = [18] . 

As we will explain below we are mainly interested in the averaged integrals of pa(z) 
and P%{z) over the target thickness W: 

f™ Piii (z)dz 

r dsc — yy ■ 

The P picture in this equation is restored by changing pa(z) — > Pi(z). In tableElthe Pd sc 
results are shown as a function of the principal and angular quantum number summed 
over the magnetic quantum number m for a 95 pm Nickel target. The differences are 
not very large, especially for the ground and lowest excited states. However, for some 
particular states the difference can be up to 20%. In figure 01 we see the discrepancy for 
the case of the 1 320) state. 

5.1. Obtaining the Break-up Probability 

Our goal is to obtain the break-up probability (P& r ) of pionium in the target. As we 
have seen in the previous sections the atoms in the target can suffer transitions between 
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Table 2. Summed P^ c = Y^m ^ds™ results in the probabilistic (P) and density 
matrix (p) pictures. The average is over W = 95 fj,m and the target material is Nickel. 



pnl 

dsc 


Pip 


1=0 


1=1 


1=2 


1=3 


1=4 


1=5 


n=l 


p 


r 


0.072854 
0.072860 












n=2 


p 


r 


0.0050676 
0.0050878 


0.008500 
0.008538 










n=3 


p 
p 


0.00087163 
0.00086909 


0.0016366 
0.0017234 


0.0020617 
0.0020250 








n=4 


p 


0.00024899 


0.0004803 


0.0006270 


0.0007326 








p 


0.00024620 


0.0005242 


0.0006445 


0.0007028 






n=5 


p 


0.000092377 


0.00018015 


0.00023838 


0.00028247 


0.00031200 






p 


0.000089072 


0.00019899 


0.00025137 


0.00027925 


0.00029343 




n=6 


p 


0.000038357 


0.000075133 


0.000099834 


0.00011906 


0.000131889 


0.00014113 




p 


0.000034640 


0.000079850 


0.000102429 


0.00011493 


0.000121041 


0.00012520 


n=7 


p 


0.000015300 


0.000029939 


0.000039634 


0.000047316 


0.000052490 


0.000057376 




p 


0.000013706 


0.000031028 


0.000039316 


0.000043479 


0.000048089 


0.000045044 



-2 

xlO 

0.35 r 




Z [\im] 

Figure 2. The solution of © and © for the |320) state. 

bound states and annihilate. However, they can also be transferred, in a collision 
with a target atom, into a continuum state. The coefficients c^j and Qik,i m accounting 
for transitions between discrete and continuum states are more difficult to compute 
than the discrete-discrete ones since the atomic form factors have a more complicate 
expression |3,. However, as shown in jHj for the case of the probabilistic picture, the 
direct calculation of break-up probability from the systems solutions is not satisfactory 
since it decreases very slowly as a function of the principal quantum number of the 
broken discrete state and only a finite number of shells (n < 7) are considered when 
solving either (JHJ) and (jEJ). We would have to guess the break up probability for any 
shell with n > 7 and make a large error in the total break-up probability determination. 
The standard strategy to obtain break-up probability consists of calculating the 
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probability of the atom to leave the target in a discrete state (Pdsc) and the probability 
of annihilation (P an h) and make use of the relation: 

1 = Pbr + Pdsc + Panh- (22) 

As both P(i sc and P an h quickly decrease with n we have an accurate result taking into 
account only those events with n < 7. A small correction will be introduced for P^ 7 - 
In the experimental conditions the atoms are not created at the target beginning 
but uniformly distributed along the target thickness. The probability that the atom 
leaves the target in a discrete state can be however linked to the solutions under ()2()j) 
initial conditions by: 

^ Pli (W - z)dz x - pij(z)dz 

p dsc-}_^ w -2^ — w — (23) 

i i 

where W stands for the target thickness (of 95 pm in our case). 

The annihilation probability is a little bit more difficult to calculate. If the atom 
is created in z , the probability that it flies to z and annihilates is given by But z 

can be any value between z and the target end W. Meanwhile, the atom is randomly 
created between and W with uniform distribution, then the annihilation probability 
is given by: 

P rw rw y r w 

Panh = Y1WJ J P^ Z ~ Z o) dzdz O = ^WJ (^-*W*)d*(24) 

Of course the probabilistic picture is restored by substituting pu(z) by P%{z) in (J23|) 
and (ED). 

As we did in (12 lj) for the P l dsc probability we can of course define the annihilation 
probability from a certain state as: 



r ; 



PLh=^ J (W-z) Pii (z)dz (25) 

where again the replacement pu(z) — > P%{z) recovers the P picture. Of course P % anh = 
for any state with k ^ 0. 

The results for the annihilation probability from the S states up to n = 7 are shown 
in table 01 and complete those of the Pd sc in table El 

In figure El we can see the dependence of P^ sc and P™ nh on the principal quantum 
number. The results have been summed over every shell bound states. We can check 
that whereas P"^ nh quickly converges to zero, and can be neglected for > 4, P dsc 
diminishes more slowly. This leads to introduce an extrapolation for P^ 7 [2/. 

+ ^ < 26 > 

where a and b are obtained by fitting P dsc at n = 5 and n — 6. The extrapolation is 
also used for n = 7 because not considering the next shells in the systems distorts this 
shell solutions. 
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Table 3. P^nh results in the P and p pictures. The average is over W — 95 /im and 
the target material is Nickel. The lifetime of pionium was assumed to be 2.9 x 10~ 15 s. 



n 


pn 

anh 


Pin 


n 


i 
i 


p 

1 








P 


U.U i zouu 


■ = 


Z 


P 


0050fi7fi 






P 


n 00^0878 
u.uuouo i o 


n = 


Q 

o 


p 








P 


u.uuuouyuy 


n = 


4 


P 


00024899 






P 


0.00024620 


n = 


5 


P 


0.000092377 






P 


0.000089072 


n = 


6 


P 


0.000038357 






P 


0.000034640 


n = 


7 


P 


0.000015300 






P 


0.000013706 




Figure 3. Dependence of P^ sc and P % anh averaged over every shell on the principal 
quantum number. The extrapolation of (|2(j|) is also shown. 

The extrapolation results are summed over n and, together with P^ 7 and P an h, 
subtracted to one to calculate the break-up probability: 

obtaining, for out particular example of 2.9 x 10~ 15 s atoms in a Ni 25 pm target: 
Pb r = 0.459254 in the probabilistic picture and Pf, r = 0.459268 in the density matrix 
formalism. The other probabilities are shown in table EJ 
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Table 4. Probability results in the P and p pictures. The average is over W = 95 
/im and the target material is Nickel. The lifetime of pionium was assumed to be 
2.9 x 1CT 15 s. 



Picture 


Pbr 


Panh 


pn<7 
dsc 


pn>7 

dsc 


P 


0.459254 


0.444536 


0.0947916 


0.001418 


P 


0.459268 


0.444575 


0.0949106 


0.001245 



6. Discussion and Conclusions 

We have checked that in the conditions of the DIRAC experiment the effect of 
the quantum interference between states does not change the result of the break-up 
probability of pionium in the target. Hence, the results obtained in the classical picture 
are accurate enough to safely perform the experimental measurement. 

The unchanged result of break-up result takes place despite the fact that for some 
discrete states, as 1 320) , the effect of interference can significantly change the population 
of the state up to 20% levels. However, the most affected states are very unpopulated 
and hence not relevant for the final results. 

The situation could change if the initial conditions were not that most atoms are 
created in the ground state. The later is non degenerated and interferences only show- 
up after a first transition. However we have checked what would happen if the initial 
conditions were that all the atoms were created in the 1 300) state and neither found a 
significant change with the probabilistic approach. A possible explanation is that while 
the interference is most likely with states with the same magnetic quantum number m, 
and comparable with the transition cross sections, the dominant transitions are those 
that increase / and m in one unit, free of interference with the father state. 
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Appendix A. The Vt matrix elements in the First Born Approximation 

Let us show how to obtain the discrete matrix elements of the Q matrix in the first Born 
approximation from the original equations of j^j. The Q operator is originally defined 
as a function of the transverse position of the atom wave functions s\ t 2- If we split the 
operator in two: 

= nW(sl,s2) + n (2) («l>^)> (a.i) 
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fi (2) (s!,£) 



r(6,si)+r*(6,s5) 



d 2 6, 



r(6,i!)r*(6, s i) d 2 b. 
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(A.2) 
(A.3) 



In the case of the ir + TT -atom the interaction operator of the Glauber theory is given 
by: 



where 



T(b, s) = 1 - exp i X {b - s/2) - z"x(6 + s/2) 



oo 



(A.4) 



(A.5) 



being £7(r) the potential of the target atoms given by the inverse Fourier transform 

of dm. 

First of all we are going to re-write f^(si, si). For that we split T(b,s) into its 
real and imaginary part: 



r(&, «i(2)) = Re r(6, si( 2 )) + i Im T(6, s 1( 



2)), 



ReT(6,s 



1(2), 



1 — COS 

1 

2 



x(6-a/2)-x(6 + s/2) 

r(6, si (2 )) r*(6, s*i(2)), 



sin 



(A.6) 

(A.7) 
(A.8) 

(A.9) 



due to the odd nature of the sin function and the even nature of x(6±s/2). Taking this 
into account we can have: 

1 



X (b-s/2)-x(b + s/2) 



Imr(6,s 1(2) ) = 

where the integral over the imaginary part goes to zero: 
Imr(6,s 1(2) )d 2 6 = 0, 



S2) 



r(b,3 1 )T*(b,3 1 )+T(b,8 2 )T*(b,3 2 ) 



d 2 b. 



Our final goal is to obtain the matrix elements ^P/L, defined as: 



n 



(1,2) 
ik, Ira 



^*(^i)^(n)^fc(^2)^m(^2)fi (1 ' 2) (sl, si) dfi dr 2 . 



In particular we can define the profile- function Tn(b): 

r«® = j i>mM7)r u (b,Z)dr, 

and its Fourier transform, the amplitude: 

2n I 



(A.10) 



(A.ll) 



(A.12) 



(A.13) 
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Tu{b) 



2ni 



e~^A u {q)d 2 q. 



It is easy to check that: 

n {2) 

ik, Ira 



~ / r«(6) Tl m {b) d 2 b=- A u (q) A* km (q) d 2 q. 
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(A.14) 

(A.15) 



To obtain an analogue of (|A.15j) for fi^ im we have to work a little bit. Of course, by 
definition: 



n 



(i) 

ik, Ira 



Jkra 




iJj*(r)Tpi(r)r(b,s)r*(b,s) df 



d z b 



+ 



Sii 




d 2 b. 



(A.16) 



r k (^ m (7)T(b,s)T*(b,s)dr 

To achieve the final result we will need the completeness equation in the form: 

5(f- ?) = V>j (r)ip* (r') , (A. 17) 

3 

which allows to express the inner integrals in (|A.16J) in terms of the profile-function 
^•(6): 

V>*(f)^(f)r(6, s) T*{b, s)df = / $(r)ih(?)6(r- r')T{b, s) T*(b, P) df df' 



E 



i(}*(r)T(b, s)ipj{r) df 



V>*(f)r*(&,sO^(f) dr' 



^T^T^b), (A. 18) 



where of course we can make the substitution: 



J2r l] (b)r lJ (b)d 2 b = J ^AMAtjifid'q, 



(A.19) 



to obtain: 



n 



(i) 

ik, Ira 



6 ~f J YlMQ)^iM d2( l+ 5 -§ J Y, A ^MAUq)d\ (A.20) 



In the Born approximation 



MO) = pU(q) 



i ( Q 
2 



(A.21) 



where we find the form factors defined in (|12p. 
Let us try to perform the sum: 



F ] 1 \- F ? -1 



Fl 



E 



W-{f) ( 



Jqr/2 _ „-iqr/2 



) i/jj(r) df 



-iqr'/2 _ e iqr'/2 



ipi{r') dr' 



V>*(r) ( f 



i< ?' r / 2 _ P -WI^\ ^p-iqr/2 _ p iqr/2 



) ^(f) df 



J V*(f) (2 - e"^ - e"^ 2 ) ^,(#0 df = 25 a - - 



(A.22) 
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From equations (f05jl . (fQ6|) . (jXTTSjt . (TP9J) . (TOTl) and (TQ2|l one can derive 
the final expressions in the Born approximation: 



O - _i_ o(2) 



(1) 

ik,lm 



n 



(2) 

ik.lm 



"k,m 



+ 
1 



2/5 2 



\U(q)\ 2 [25 hl -Fl(q)-Fl(-q)] d 2 q + 
J \U(q)\ 2 [24, m - F™{q) - J^(-g)] d 2 g, 

M ( € 



P(q)\ 



d 2 q. 



(A.23) 



(ED 



(Hi 
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